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ABSTRACT 

The  development   of  a  numerical   1000-mb   prognosis  by   means 
of   a   500   to    1000-mb   thickness  forecast  is   investigated.      The 
thickness   forecast  is  based,    in  part,  on  the  500-mb  barotropic 
prognosis  as   developed  by   Fleet    Numerical  Weather   Facility, 
Monterey. 

Forecasts   were  prepared  for   three  days   during  May   1966, 
•with  various   amounts   of   smoothing  of   the  component  fields. 
The   model   showed  promise  of  producing  a   rapid,    usable  prognosis; 
however  further   testing  is   needed  with  the  probable  addition  of 
some  climatic   constraints.      It  is   also  hoped   that   the  yet  untried 
diabatic  heating  and  surface  friction  terms   will  provide   some 
further  improvement. 
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SUBSCRIPTS 

a  air  at  surface 
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1 .      Introduction 

With  the  development  of  the  high  speed  digital  computer,  the 
problem  of  weather  prognosis  has  evolved  from  one  of  qualitative 
reasoning  by  an  experienced  forecaster  to  one  of  objective 
computations  based  on  mathematical  equations  representing  the 
motions  of  the  atmosphere.    At  this  time,  daily  numerical  analysis 
and  prognostic  charts  are  being  produced  at  Fleet  Numerical  Weather 
Facility,  Monterey  (FNWF) ,  and  rapidly  disseminated  to  fleet  units 
around  the  world.    The  entire  field  of  meteorology  is  in  a  rapid  state 
of  change  with  more  and  more  products  being  produced  numerically. 

The  most  successful  numerical  prognosis  has  been  at  the  500  mb 
level.    If  a  successful  numerical  isobaric  thickness  forecast  could 
be  developed,  the  success  of  the  500-mb  forecast  could  then  be 
reflected  to  higher  and  lower  levels.    The  objective  of  this  research 
is  to  obtain  a  dependable  24-hour  numerical  forecast  of  thickness 
for  the  layer  between  500-mb  and  1000-mb,  and  then  using  FNWF's 
500-mb  barotropic  forecast  obtain  a  1000-mb  height  prognosis. 

The  thickness  of  an  air  layer  at  any  given  time  and  place  is 
dependent  upon  the  mean  temperature  of  the  air  column.    If  we  look 
at  thickness  change  as  a  result  of  the  advection  of  thickness  lines 
only,  it  appears  that  there  is  no    heat  loss  or  gain  as  the  air  moves 
across  the  surface  of  the  earth  nor  are  there  vertical  motions.    This 
is  not  true,  of  course,   since  the  air  gains  or  loses  heat  through 
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radiative,  dynamical,  and  conductive  processes,  and  through  latent 
heat  exchange  resulting  from  phase  changes  of  the  water  constituent. 

In  the  model  used  in  this  study,  in  addition  to  advection,  an 
empirical  heating  term  is  developed  to  acdount  for  heating  or  cooling 
by  evaporation  or  condensation    and  for  turbulent  transfer  of 
sensible  heat  over  ocean  surfaces. 

To  account  for  the  dynamic  processes,  vertical  velocity  is 
assumed  to  be  a  second  order  function  of  pressure.    The  1000-mb 
vertical  velocity  is  composed  of  a  terrain  upslope  term  and  a 
friction  term.    Heat  gained  or  lost  by  radiative  processes  was 
ignored  as  neglibible  over  the  24-hour  period. 
2.      Development  of  the  Prognostic  Model 

A  two  layer  model  is  used  in  the  development  of  the  prognostic 
equations.    The  lower  troposphere  is  divided  into  two  layers  as 
illustrated  in  Fig.   1. 

By  use  of  the  equation  of  state  and  the  hydrostatic  equation, 
the  thermodynamic  equation  may  be  reduced  to  the  following  form 

Replacing  the  vertical  derivatives  above  by  finite  difference  and 
simplifying  leads  to 
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Where    A£=h    ,  the  thickness. 

In  the  layer  1000  to  500  mb  it  has  been  noted  that  6)  can 
generally  be  approximated  by  a  second  order  equation  in  pressure 

^  =  aowt)  feeds)*  -+  b(v,t)(^t>)  +  io,o         (3) 

By  means  of  equation(2)it  may  be  shown  that 


K=  %*-  +  3/4b  +wio 


(4) 


6Mb  r  ^ol  -\.  j/^b  +WW  (5) 

Taking  the  vorticity  equation  at  500  mb  in  the  form 
and  substituting  from  equation   (6)  gives 
Equation  (2)  is  now  evaluated  for  layer  2  giving 

aKv+V.\7k    4.(a%(-S»-,,U  -(birULMlQ 

The  above  equation  may  be'  -written 


(7) 
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Where 


*■      \Gff6n-  Af       / 

Combing  equations  (5)  and  (8)  and  simplifying  leads  to 

With  the  geostrophic  assumption  for  relative  vorticity 

i  -  wan 

it  follows  that 

TO-vn^.Ki      (id 

Subtracting  equation  (ll)  from  (10)  yields 


However, 


3  /P  /,  ^>*>  X 


P,  »    -fe  '      AA*  =  ft 

thus  equation  (12)  becomes 

W^h*-  =V,*VVV*^4l^\i  -%(a±b)    (13) 
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Substituting  equations  (7)  and  (9)  into  equation  (13)  and  expressing 
advective  terms  in  Jacobian  form  we  obtain  the  following: 

where  ,    ''"  ' 

and  ' 

The  prognostic  thickness  equation  is  thus  determined  far  the  lower 
layer.    This  equation  is  of  the  Helmholtz  type,  the  solution  of 
which  may  readily  be  obtained  by  use  of  a  digital  computer. 

We 'Will  now  obtain  a  prognostic  equation  for  the  thickness 
of  the  upper  layer. 

Evaluating  equation  (2)  for  the  upper  level  centered  at  625  mb 
results  in  the  following: 

where  ^_  x 

Vcft,'         ^T  ' 

Substituting  from  equations  (4^  (7)  and  (9  )  and  simplifying,  we 
obtain 

<*•  -iy-VK,-^^-!^^)-^.^)    (16) 
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Substitution  of  equation  (16)  into  equation  (15),  and  expressing 
in  Jacobian  form  provides  the  thickness  tendency  (prognostic) 
equation  for  the  upper  layer  of  the  model. 


^  : 


where 

V-Vh,  =  %/4  J(P-,;Dr) 

III.    Development  of  the  heating  term. 

Initially  it  will  be  assumed  that  the  heating  rate  is  expressed 
by  an  equation  of  the  following  format 


6-A|(V)tTs-TJ^j 


(18) 


Over  land  the  heating  term  is  quite  small,  and  we  shall 
consider  it  as  negligible.     Over  water,  excluding  radiation,  the 
main  contributions  to  the  heating  term  are  due  to  sensible  and 
latent  heat  exchange  between  the  sea  and  the  atmosphere.    It  is 
of  course  realized  that  latent  heat  gained  through  evaporation  may 
not  be  utilized  by  the  atmosphere  at  the  exact  instant  and  position 
that  it  is  obtained.     However,  in  the  relationship  derived,  an 
empirical  adjustment  constant  is  available  should  the  values 
calculated  be  determined  to  be  somewhat  too  small  or  large. 
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According  to  Laevastu   C  3  J    the  heat  flux  from  the  ocean  to 
the  atmosphere  due  to  sensible  heat  transfer  and  evaporation  or 
condensation  may  be  expressed  as  follows: 

Qw  -  g  (M(. +crn  vf)(Ts  -T. )  ( £££)      ,  19 , 

Qfe    *    L.(0.au*0/77\£)(.18es-eo^  (20) 

£4 


(21) 


QFC  =  .077W(es-(fOL 

a4  (22) 

for  CTs-To.)<0 


It  is  apparent  that  the  contribution  from  the  latter  equation 
is  small  and  hence  will  be  considered  zero. 

By  equating  the  relationship  (19)  given  by  Laevastu  for 
sensible  heat  flux  with  the  integral  of  equation  (18)  with  respect 
to  mass  per  unit  area,  we  may  determine  a  reasonable  value  for 
the  constant  "A" . 

Integrating  (18)  with  respect  to  mass  per  unit  area  gives 

2 


qf  =  ]hmvs-x)^fi 


where 

Introducing  the  hydrostatic  relationship,  p  d?:-  -d/p  and  integrating 

r  5 

with  respect  to  pressure,  and  then  equating  this  term  to  the  heat 

flux  (19)  yields 
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A4(VsX*o-4»ft(Ts-TJ  -  5l(a2b+o.-?7vO(Ts-TJ(ija^ 

3<*  W ~     ^  ^ 

Solving  this  equation  for  A,  and  placing  ,/(U)-  {0>lU>+b.1M) 

leads  to  a  value  of  63.5  for  A.    The  resulting  term  for  the  sensible 

heat  contribution  becomes 

QH  =   1&.S  \&.V*  +ai1V0(Ts-Tfc)^j*j*     (^J  <23> 

Furthermore ,  since  the  wind  used  in  the  prognostic  model  will 
be  the  geostrophic  wind  at  1000  mb  which  is  about  130  to  140  per- 
cent greater  than  the  surface  wind,  the  constant  0.77  in  equation  (23) 
was  reduced  to  0.5.    This  constant  may  be  refined  to  further  enhance 
the  accuracy  of  the  forecast. 

The  final  form  for  this  term  thus  becomes 


QH  =  (obS  (p^+osv^KTs-rj^f  (24) 

for     (Ts-TO^O 


Bowen  (1926)  related  the  convective  transfer  of  sensible  heat 
to  evaporation  by  the  following  relationship 


0e  (es  -eO  1000 

A.   H.   Gordon    £2]     (1952)  using  meteorological  data  from 
British  ocean  weather  ships  in  the  North  Atlantic  related  the  mean 


(25) 
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values  of  vapor  pressure  differences  (e   -e  )  to  individual  whole 

o  ci 

degree  temperature  differences  of  (Ts-Ta)  for  wind  forces  0-3,  4-6 
and  greater  than  6. 

Roll    C4j     using  this  data  computed  the  Bowen  ratio  from 
corresponding  average  temperatures  and  water  vapor  differences 
(Ta-Ts),   (es-ea)  for  wind  forces  Beaufort    4    and    8.     Only  a  very 
slight  increase  in  the  Bowen  ratio  was  observed  for  force  8  winds. 
The  results  of  those  computations  are  presented  in  Fig.  2. 

The  curve  of  Fig.  2  was  statistically  fitted  using  a  linear 
rational  transformation  and  solving  for  the  Bowen  ratio  in  terms  of 
temperature  difference  (Ts-Ta)  with  the  following  result: 

P>  =    Q-*      =     *  3ATs-T*-)  (26) 

As  was  previously  indicated,  the  diabatic  heating  term  is  to 
be  comprised  only  of  the  two  elements    QH       and  QE    ,  hence  the 
use  of  relationship  £6)  enables  us  to  eliminate  the  vapor  pressure 
differences  as  a  parameter  in  the  forecast.    Thus   Q    (Total  Diabatic 
Heating  Rate)  =  Q^.  +  Qe 

The  total  heating  rate  determined  from  equations  (24)  and  (27) 
is  thus 
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From  equations  (2)  and  (18) ,  the  rate  of  change  of  thickness 
due  to  diabatic  heating  in  the  lower  and  upper  layers ,  which  will 
be  termed  Q2  and  Q-^  respectively,  may  be  expressed  as  follows 

-  y»f/Aiiv)(T.-To  (29) 

Q,=  ^^M(Vs)(T,-TJ  (30) 

or         Q.  =     5H  Q^  (3D 


\ 

We  have  thus  established  a  simplified  distribution  of  the 
effects  of  the  diabatic  heating. 

From  equations  (28)  and  (31)  the  rate  of  change  of  thickness 
due  to  diabatic  heating  in  the  lower  layer  is  then 

=    \S.±(\  +  .l3.(Ti~T^)){.^+.S\/p,o)    (g^)         (32) 

£or    (Ts-TJ?o 

Similarly,  it  may  be  shown  that 

Q^=0.Xl(vVlo)(Ts-To.)  (33) 

for        (Ts-TO<  O 
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In  comparison  with  equation  (32) ,  it  is  apparent  that  the 
magnitude  of  the  thickness  change  due  to  this  term  will  be  quite 
small.     However,  the  constant  .27  may  be  increased  if  this 
contribution  is  deemed  too  small. 
IV.    Development  of  the  terrain  pressure  and  friction  term 

According  to  Berkof  ski  and  Bertoni    [  ll     ,  the  vertical  velocity- 
produced  by  terrain  may  be  expressed  as  follows: 


ut  =  -f^VT'V£T 


(34) 


Making  the  geostrophic  assumption,  and  also  assuming  that 
the  1000-mb  wind  is  valid  for  all  levels  of  terrain,  this  relationship 
may  be  expressed  as  follows: 

w--e0i\j(D(o,*T)  (35) 

It  is  apparent  that  this  term  may  be  somewhat  small  in  cases 
of  high  terrain  due  to  the  fact  that  the  geostrophic  wind  at  1000  mb 
is  used  rather  than  an  integrated  wind  over  the  entire  height  of  the 
terrain.    However,  in  nature,  the  wind  will  normally  tend  partially 
to  parallel  the  terrain  contours ,  hence  this  reduction  will  somewhat 
compensate  for  this  deficit. 

It  was  originally  intended  to  use  a  more  refined  function  for 
terrain  induced  vertical  velocities  and  to  include  the  effects  of 
friction.    Time  limitations  however  precluded  any  further  work  on 
this  term. 
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It  is  felt  that  further  expansion  of  this  term  or  the  inclusion 

of  the  terrain  term  now  calculated  by  FNWF  into  the  program  will 

render  some  improvement  to  the  prognosis. 

V.      Programming  the  Prognostic  Equations  for  the  Control 
Data  1604  Computer 

The  prognostic  equations  were  put  into  finite-difference  form 
and  programmed  for  the  Control  Data  1604  computer  used  in 
conjunction  with  the  FE880  Univac  storage  drum.    The  1604  computer 
is  a  stored  program,  general  purpose  digital  computer  with  a  storage 
capacity  of  32,768  words  in  two  indpendent  storage  units.    The 
FE880  storage  drum  is  connected  to  the  1604  computer  by  a  high 
speed  channel  which  conveys  information  at  a  rate  of  3,000,000 
bits  per  second.    The  drum  capacity  is  256  blocks  of  5,120 
48  bit  words  per  block. 

The  program  was  written  in  symbolic  coded  relocatable 
assembly  program  (SCRAP)  for  a  63  by  63  hemispheric  grid,  with 
one  hourly  time  steps . 

An  option  was  written  into  the  program  to  allow  for  equating 
the  vertical  velocity  due  to  terrain  and  the  heating  term  to  zero 
in  order  to  ascertain  their  respective  effects  on  the  prognosis. 

The  750-mb  layer  is  not  considered  a  standard  layer  in  the 
atmosphere,   so  no  synoptic  information  is  available  for  the 
initialization  of  the  program .    To  account  for  this ,  the  following 
equation  developed  from  the  thickness  equation  was  derived. 
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Diso  =  Dioo  ■*  Kern  ~  If8"73  (T«»-  Ti^  +  3oiT' -3        1 36 » 

where  h  t77,-  i-s  tne  standard  thickness  of  the  layers  700  to  750  mb. 

Once  the  initialization  of  the  program  is  completed  the      750  mb 
height  may  be  computed  using  the  forecast  thickness  for  the  upper 
level  and  FNWF's  500-mb  barotropic  prognosis. 

The  absolute  vorticities  (^]s  ^  /ft    )  are  computed  with  FNWF's 
subroutine  SAR.    This  routine  computes  the  absolute  vorticity  over 
a  63  by  63  grid  using  the  following  relationships 

\  =  4,  Sine  <0    -r-  a^  (A.+A^+As  -^A^-^Aq) 

I  Sine    6{ 

Where  al  =  .0746711 

a2  =  .2704348 

-17 

To  use  this  routine,  all  data  must  be  scaled  2        .    The 

absolute  vorticity  is  equated  to  zero  for  all  points  south  of  10  N 
latitude . 

After  computing  the  absolute  vorticities  and  reading  in  the 
various  D  fields,  the  Jacobians  may  be  computed  utilizing  Fleet 
Numerical  Weather  Facility  standard  subroutine  SAB.     The 
subroutine  computes  the  Jacobian  in  the  following  manner 

A  =  T(e>,0 

=  (A„(3)(AyC)-(A^6)(AxC) 

The  program  inserts  zeroes  in  field  A  at  the  boundaries. 
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The  next  term  computed  is  the  V  L^V^     .    This  term  is 
computed  by  Fleet  Numerical  Weather  Facility1  s  subroutine  SAD , 
which  computes  as  follows: 


=  6t 


+1)i^H,iAri+^M"^i 


This  program  also  zeroes  the  field  at  the  boundaries. 

To  compute  the  heating  term,  Qn  we  must  obtain  an  hourly 
surface  temperature  based  on  the  forecast  thickness  and  also 
compute  a  1000-mb  geostrophic  wind.    The  sea  surface  temperature 
was  assumed  to  remain  constant  throughout  the  period  of  the 
forecast. 

The  equation  developed  to  determine  an  hourly  surface 
temperature  as  derived  from  the  predicted  thickness  fields  is  as 
follows: 

X-c  =3.4-lS"x  I0"4[(i0»  hK-133h.>jU,  (8-'»"ftPJ+a.<7lj|  (37) 

After  Q9  is  calculated,  each  point  in  the  grid  is  checked  to 
determine  whether  it  is  a  land  or  sea  point.    At  land  points  Q9 
is  equated  to  zero;  at  sea  points,  as  calculated. 

The  terrain  term  W-,  q  is  computed  using  the  Jacobian  subroutine 
and  also  checked  against  the  land/sea  table.    In  this  case,   sea 
points  are  zeroed  and  land  points  left  unchanged. 


24 


The  terms  now  calculated  are  collected  and  make  up  the  forcing 
function  for  the  Helmholtz  equation  {Ik)-    This  equation  is  solved 
using  FNWF's  subroutine  SAE,  providing    Z\  rL  for  the  first  hour  of 
the  prognosis. 

To  solve  the  Helmholtz  equation  ,\rA'BA  s  C  ,  the  program  scans 
through  field  A  laterally  from  the  lower  left  point  to  the  upper  right 
point  using  the  extrapolated  Liebmann  method  of  computing  a  new 
point  value.    The  boundary  conditions  have  already  been  determined 
and  stored  in  the  A  field,     the  guess  field,which  is  originally  set 
to  zero. 

A*!'    =  A', :   +  A  f  ^  -  (AB):,i  -CoiT 
1/1  LA     ¥  L       a-*  Btf+i J 

Ka  =  ^  P^U  -\A6^-CtltV 

When  for  all  points    I  Kv.,^1^-  €     where  £is  the  criterion  for 
convergence  or  the  point  where  the  residual  R  is  considered 
liquidated,  the  field  is  assumed  relaxed  and  a  proper  solution  has 
been  obtained  at  each  grid  point. 

When  A  ho  has  been  determined,  all  the  parameters  required 
for  the  computation  of  Ahj  are  available  for  use  in  equation  (17) . 

The  change  in  height  at  1000  mb  is  then  determined  by 
subtraction  of  A  h^  plus  A  h2  from  the  500-mb  height  change  as 
computed  by  the  FNWF  barotropic  model. 


25 


6.      Results  and   Conclusions 

Twenty  four  hour   thickness   forecasts   were  produced  from 
initial  data   on  20,    23   and  29  May  1966.      The  initial  program 
contained  a  number   of   smoothing  passes  at   many  stages   of   the 
prediction  procedure.      As  a  result  the  initial  thickness  prognosis 
■was   very   smooth,    almost   zonal,    as   shown  in  Figure  4-      Compar- 
ison of   this   prognosis   with   the  initial  and  verifying   maps, 
Figures  3   &nd  6»    clearly   shows   that  oversmoothing  was  a   major 
source   of   error.      In  figures   5f    10   and   13   for   the  20th,    23rd 
and  29th   of  May  the  results   were  slightly   more   encouraging. 
Here  smoothing  has  been  reduced  somewhat  and  lower   limits 
were  placed   on  absolute  vorticity  and  stability.      Still  further 
elimination  of   smoothing  throughout  the  prediction  procedure 
gave   the  thickness  and   1000-mb   prognosis  for   the  29th  in 
Figures    15   and   19,    which  show   genuine  promise  as   a   useable 
forecast.      No    1000-mb  analysis    existed  for   the  21st  and  24th. 

Generally,    the  pressure   systems   moved  realistically  and 
the  final  positions   were  good.      However,    some   excess   move- 
ment was   observed  which   might  be  partially  attributed  to   the 
smoothing   still  remaining  in  the   model.      Some   overdevelopement 
of  both  high  and  low  pressure  systems   in  Figure   15   was   also 
observed.      However   FNWF  has   a  filtering   subroutine  which  in 
addition  to  filtering   out  the  ultra   short   waves   reduces   the 
amplitude   of   these   centers.      This   subroutine  diminishes   the 
effect   of   simply   spreading   the  patterns   to  adjacent  grid  points 
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which  is   inherent   in   smoothing  subroutines.       Smoothing  can  also 
give    rise   to   fictitious    movement   of         pressure   systems. 

By  replacing   the   smoothing   with  a   filtering   pass   some  of 
the  difficulties   of   excess   movement  and  overdevelopement   could 
be   overcome. 

Another  problem   is   the   tendency  of   some   thermal   systems, 
known  for   their   persistance,    to    move   with  the  general  circulation, 
an  example  is   the  heat  low   in  Southwestern  United   States. 
Periodic   comparison   with   climatology  could  be  fruitful  in  reducing 
this   tendency. 

Due  to  a  lack  of   available  computer   time  little   consideration 
could  be  given  to   the  diabatic  heating  and   terrain  terms.      Further 
refinement  and  inclusion   of   these   terms   could  result   in  some 
improvement   of   the  prognosis. 

The  program   running  time  on  the   1604   computer  for  a  24- 
hour  forecast  was   about  30    minutes,    including  hourly   writing 
on   magnetic   tape   of  a  number   of  fields   which  would  be  unnecessary 
in  an  operational   model. 

The   FNWF   thickness  forecast  in  Figure   16,    which  is  based 
on  simple  horizontal  advection  and  numerous   constraints,    is 
obviously      superior   to   the  prognosis   resulting  from    the   model 
under  investigation...    However,    the  superiority  is   a  result   of 
several  years   of   refinement  and   empirical  adjustment.      Time 
did  not   permit   such   thorough   investigation  of   this   research   model. 

In  conclusion,    it  has  been  determined   that   the  assumptions 

27 


used  in  the   construction   of    this    model  are  reasonable.      With 
proper  adjustment  and  refinement   this    model   should  yield  a  good 
forecast   of   thickness.      This,    in  turn,    should  produce  a  good 
surface  prognosis. 
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Figure   1.      Two  Layer   division   of   the  atmosphere  between 
1000-mb  and   500-mb.      h^   and  hg   are   the  thickness  of   the 
layers   p^    to   p^   and  p10    to   p«   respectively. 
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Figure  2.       Bowen   Ratio   vs.    Sea-Air   Temperature    Difference 
(dashed   curve) 

Linear   rational   transformation  for   positive  values, 
Equation    (26)    (Solid   curve). 
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